{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The point and line-segment plotting provided by Datashader can be put together in different ways to visualize specific types of data. For instance, network graph data, i.e., networks of nodes connected by edges, can very naturally be represented by points and lines.  Here we will show examples of using Datashader's graph-specific plotting tools, focusing on how to visualize very large graphs while allowing any portion of the rendering pipeline to replaced with components suitable for specific problems.\n",
    "\n",
    "First, we'll import the packages we are using and demonstrating here.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import math\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "import datashader as ds\n",
    "import datashader.transfer_functions as tf\n",
    "from datashader.layout import random_layout, circular_layout, forceatlas2_layout\n",
    "from datashader.bundling import connect_edges, hammer_bundle\n",
    "\n",
    "from itertools import chain"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Graph (node) layout\n",
    "\n",
    "Some graph data is inherently spatial, such as connections between geographic locations, and these graphs can simply be plotted by connecting each location with line segments. However, most graphs are more abstract, with nodes having no natural position in space, and so they require a \"layout\" operation to choose a 2D location for each node before the graph can be visualized.  Unfortunately, choosing such locations is an [open-ended problem involving a complex set of tradeoffs and complications](http://www.hiveplot.com).\n",
    "\n",
    "Datashader provides a few tools for doing graph layout, while also working with external layout tools. As a first example, let's generate a random graph, with 100 points normally distributed around the origin and 20000 random connections between them:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.random.seed(0)\n",
    "n=100\n",
    "m=20000\n",
    "\n",
    "nodes = pd.DataFrame([\"node\"+str(i) for i in range(n)], columns=['name'])\n",
    "nodes.tail()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "edges = pd.DataFrame(np.random.randint(0,len(nodes), size=(m, 2)),\n",
    "                     columns=['source', 'target'])\n",
    "edges.tail()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here you can see that the nodes list is a columnar dataframe with an index value and name for every node.  The edges list is a columnar dataframe listing the index of the source and target in the nodes dataframe.  \n",
    "\n",
    "To make this abstract graph plottable, we'll need to choose an x,y location for each node. There are two simple and fast layout algorithms included:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "circular  = circular_layout(nodes, uniform=False)\n",
    "randomloc = random_layout(nodes)\n",
    "randomloc.tail()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cvsopts = dict(plot_height=400, plot_width=400)\n",
    "\n",
    "def nodesplot(nodes, name=None, canvas=None, cat=None):\n",
    "    canvas = ds.Canvas(**cvsopts) if canvas is None else canvas\n",
    "    aggregator=None if cat is None else ds.count_cat(cat)\n",
    "    agg=canvas.points(nodes,'x','y',aggregator)\n",
    "    return tf.spread(tf.shade(agg, cmap=[\"#FF3333\"]), px=3, name=name)\n",
    "\n",
    "tf.Images(nodesplot(randomloc,\"Random layout\"),\n",
    "          nodesplot(circular, \"Circular layout\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The circular layout provides an option to distribute the nodes randomly along the circle or evenly, and here we've chosen the former.\n",
    "\n",
    "The two layouts above ignore the connectivity structure of the graph, focusing only on the nodes. The [ForceAtlas2](http://journals.plos.org/plosone/article/file?id=10.1371/journal.pone.0098679&type=printable) algorithm is a more complex approach that treats connections like physical forces (a force-directed approach) in order to construct a layout for the nodes based on the network connectivity:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%time forcedirected = forceatlas2_layout(nodes, edges)\n",
    "tf.Images(nodesplot(forcedirected, \"ForceAtlas2 layout\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This algorithm is designed to place densely connected nodes closer to each other, but of course we will only be able to evaluate how well it has done so once we plot edges (below).\n",
    "\n",
    "## Edge rendering/bundling\n",
    "\n",
    "Assuming that we have a suitable layout for the nodes, we can now plot the connections between them.  There are currently two bundling algorithms provided: drawing a line directly between any connected nodes (``connect_edges``), and an iterative \"bundling\" algorithm ``hammer_bundle`` (a variant of [Hurter, Ersoy, & Telea, ECV-2012](http://www.cs.rug.nl/~alext/PAPERS/EuroVis12/kdeeb.pdf)) that allows edges to curve and then groups nearby ones together to help convey structure. Rendering direct connections should be very quick, even for large graphs, but bundling can be quite computationally intensive."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def edgesplot(edges, name=None, canvas=None):\n",
    "    canvas = ds.Canvas(**cvsopts) if canvas is None else canvas\n",
    "    return tf.shade(canvas.line(edges, 'x','y', agg=ds.count()), name=name)\n",
    "    \n",
    "def graphplot(nodes, edges, name=\"\", canvas=None, cat=None):\n",
    "    if canvas is None:\n",
    "        xr = nodes.x.min(), nodes.x.max()\n",
    "        yr = nodes.y.min(), nodes.y.max()\n",
    "        canvas = ds.Canvas(x_range=xr, y_range=yr, **cvsopts)\n",
    "        \n",
    "    np = nodesplot(nodes, name + \" nodes\", canvas, cat)\n",
    "    ep = edgesplot(edges, name + \" edges\", canvas)\n",
    "    return tf.stack(ep, np, how=\"over\", name=name)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cd = circular\n",
    "fd = forcedirected\n",
    "\n",
    "%time cd_d = graphplot(cd, connect_edges(cd,edges), \"Circular layout\")\n",
    "%time fd_d = graphplot(fd, connect_edges(fd,edges), \"Force-directed\") \n",
    "%time cd_b = graphplot(cd, hammer_bundle(cd,edges), \"Circular layout, bundled\")\n",
    "%time fd_b = graphplot(fd, hammer_bundle(fd,edges), \"Force-directed, bundled\") \n",
    "\n",
    "tf.Images(cd_d,fd_d,cd_b,fd_b).cols(2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The four examples above plot the same network structure by either connecting the nodes directly with lines or bundling the connections, and by using a random layout or a force-directed layout.  As you can see, these options have a big effect on the resulting visualization. \n",
    "\n",
    "Here we'll look more closely at the bundling algorithm, using a simple example where we know the structure: a single node at the center, with random points on a circle around it that connect to the central node (a star graph topology):\n",
    "\n",
    "<!-- def circle(r,n): return [(math.cos(2*math.pi/n*x)*r,math.sin(2*math.pi/n*x)*r) for x in range(0,n)] -->"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 75\n",
    "np.random.seed(0)\n",
    "x = np.random.random(n)\n",
    "\n",
    "snodes = pd.DataFrame(np.stack((np.cos(2*math.pi*x),\n",
    "                               np.sin(2*math.pi*x))).T, columns=['x','y'])\n",
    "snodes.iloc[0] = (0.0,0.0)\n",
    "sedges = pd.DataFrame(list(zip((range(1,n)),[0]*n)),columns=['source', 'target'])\n",
    "star = snodes,sedges"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tf.Images(graphplot(snodes, connect_edges(*star),\"Star\"),\n",
    "          graphplot(snodes, hammer_bundle(*star),\"Star bundled\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here you can see the bundling algorithm forms groups of nearby connnections, which helps make the structure at a particular scale clear.  The scale of this structure, i.e., how much bundling is done, is determined by an effective \"bandwidth\", which is a combination of an `initial_bandwidth` parameter and a `decay` time constant for annealing this bandwidth over time:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "grid =  [graphplot(snodes, \n",
    "                   hammer_bundle(*star, iterations=5, decay=decay, initial_bandwidth=bw),\n",
    "                                 \"d={:0.2f}, bw={:0.2f}\".format(decay, bw))\n",
    "    for decay in [0.1, 0.25, 0.5, 0.9] for bw    in [0.1, 0.2, 0.5, 1]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tf.Images(*grid).cols(4)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Starting from the bottom left and moving diagonally to the upper right, the scale of the bundling increases along a diagonal to the upper right, with higher initial bandwidth and higher decay time constants leading to larger-scale bundling.  For the largest decay time constant, the algorithm has failed to converge for large initial bandwidths (the bw 0.5 and 1.0 plots on the bottom row), because the algorithm stops at a specified maximum `iterations`, rather than reaching a fully organized state.\n",
    "\n",
    "Of course, even when the algorithm does converge, larger amounts of bundling can magnify small amounts of clumping over large scales, which may or may not be relevant to the questions being asked of this data, so it is important to set these parameters appropriately for the types of structures of interest.\n",
    "\n",
    "<!--\n",
    "max_iterations=10\n",
    "hmap = hv.HoloMap({(it, bw, decay): hv.Curve(hammer_bundle(nodes.data, edges.data,\n",
    "                         decay=decay, initial_bandwidth=bw, iterations=it))\n",
    "                     for decay in [0.1, 0.25, 0.5, 1, 2] \n",
    "                     for bw in [0.1, 0.2, 0.5, 1] \n",
    "                     for it in range(max_iterations)},\n",
    "                   kdims=['Iteration', 'Initial bandwidth', 'Decay'])\n",
    "    \n",
    "nodes_ds = datashade(nodes,cmap=[\"cyan\"])\n",
    "datashade(hmap.grid(['Initial bandwidth', 'Decay']), **sz).map(lambda e_ds: e_ds * nodes, hv.DynamicMap)\n",
    "-->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Graphs with categories\n",
    "\n",
    "One of the main uses for visualizations of large graphs is to examine the connectivity patterns from nodes of different categories. Let's consider an artificial example with four groups of highly interconnected nodes:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.random.seed(1)\n",
    "cats,n,m = 4,80,1000\n",
    "\n",
    "cnodes = pd.concat([\n",
    "           pd.DataFrame.from_records([(\"node\"+str(i+100*c),\"c\"+str(c)) for i in range(n)], \n",
    "                        columns=['name','cat']) \n",
    "             for c in range(cats)], ignore_index=True)\n",
    "cnodes.cat=cnodes.cat.astype('category')\n",
    "\n",
    "cedges = pd.concat([\n",
    "           pd.DataFrame(np.random.randint(n*c,n*(c+1), size=(m, 2)), \n",
    "                        columns=['source', 'target'])\n",
    "         for c in range(cats)], ignore_index=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The ``cnodes`` and ``cedges`` data structures form a graph that has clear structure not visible in a random layout, but is easily extracted using the force-directed approach:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rd = random_layout(     cnodes, cedges)\n",
    "fd = forceatlas2_layout(cnodes, cedges)\n",
    "\n",
    "%time rd_d = graphplot(rd, connect_edges(rd,cedges), \"Random layout\",          cat=\"cat\")\n",
    "%time fd_d = graphplot(fd, connect_edges(fd,cedges), \"Force-directed\",         cat=\"cat\") \n",
    "%time rd_b = graphplot(rd, hammer_bundle(rd,cedges), \"Random layout, bundled\", cat=\"cat\")\n",
    "%time fd_b = graphplot(fd, hammer_bundle(fd,cedges), \"Force-directed, bundled\",cat=\"cat\") \n",
    "\n",
    "tf.Images(rd_d,fd_d,rd_b,fd_b).cols(2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As you can see, the highly interconnected subgroups are laid out in separate locations in the plane, mostly non-overlapping, allowing these groups to be detected visually in a way that they aren't in a random layout, with or without bundling."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using graphs from NetworkX\n",
    "\n",
    "The above examples constructed networks by hand.  A convenient way to get access to a large number of [graph types](https://networkx.github.io/documentation/stable/reference/generators.html) is the separate [NetworkX](https://networkx.readthedocs.io) package.  Here, we will select several standard graph structures, lay them each out in the same fixed circular shape using NetworkX, and then show how they will appear without bundling, with moderate levels of bundling, and with high amounts of bundling."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import networkx as nx\n",
    "\n",
    "def ng(graph,name):\n",
    "    graph.name = name\n",
    "    return graph\n",
    "\n",
    "def nx_layout(graph):\n",
    "    layout = nx.circular_layout(graph)\n",
    "    data = [[node]+layout[node].tolist() for node in graph.nodes]\n",
    "\n",
    "    nodes = pd.DataFrame(data, columns=['id', 'x', 'y'])\n",
    "    nodes.set_index('id', inplace=True)\n",
    "\n",
    "    edges = pd.DataFrame(list(graph.edges), columns=['source', 'target'])\n",
    "    return nodes, edges\n",
    "\n",
    "def nx_plot(graph, name=\"\"):\n",
    "    print(graph.name, len(graph.edges))\n",
    "    nodes, edges = nx_layout(graph)\n",
    "    \n",
    "    direct = connect_edges(nodes, edges)\n",
    "    bundled_bw005 = hammer_bundle(nodes, edges)\n",
    "    bundled_bw030 = hammer_bundle(nodes, edges, initial_bandwidth=0.30)\n",
    "\n",
    "    return [graphplot(nodes, direct,         graph.name),\n",
    "            graphplot(nodes, bundled_bw005, \"Bundled bw=0.05\"),\n",
    "            graphplot(nodes, bundled_bw030, \"Bundled bw=0.30\")]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "n=50\n",
    "plots = [nx_plot(g) for g in\n",
    "           [ng(nx.complete_graph(n),        name=\"Complete\"), \n",
    "            ng(nx.lollipop_graph(n, 5),     name=\"Lollipop\"),     \n",
    "            ng(nx.barbell_graph(n,2),       name=\"Barbell\"),\n",
    "            ng(nx.ladder_graph(n),          name=\"Ladder\"),   \n",
    "            ng(nx.circular_ladder_graph(n), name=\"Circular Ladder\"), \n",
    "            ng(nx.star_graph(n),            name=\"Star\"),\n",
    "            ng(nx.cycle_graph(n),           name=\"Cycle\")]]\n",
    "\n",
    "tf.Images(*chain.from_iterable(plots)).cols(3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As you can see, both bundled and unbundled representations reflect important aspects of the graph structure, but the bundling results do depend on the parameters chosen.  Bundling is also very computationally expensive; nearly all of the time taken to render these plots is for the bundling step.\n",
    "\n",
    "Note that the `star_graph` example above differs from the one in the previous sections, in that all nodes here connect to a node on the outer circle instead of one in the center, which shows clearly how the layout can affect the resulting visualization."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Interactive graphs\n",
    "\n",
    "The above plots all show static images of nodes and edges, with optional category information, but there's no way to see the specific identity of individual nodes.  With small numbers of nodes you can try coloring them to convey identity, but in general the only practical way to reveal identity of nodes or edges is typically interactively, as a user inspects individual items.  Thus interactive plots are often necessary for doing any exploration of real-world graph data.\n",
    "\n",
    "The simplest way to work with interactive datashaded graphs is to use [HoloViews](http://holoviews.org), which includes specific support for [plotting graphs with and without Datashader](http://holoviews.org/user_guide/Network_Graphs.html):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import holoviews.operation.datashader as hd\n",
    "import holoviews as hv\n",
    "hv.extension(\"bokeh\")\n",
    "\n",
    "circle = hv.Graph(edges, label='Bokeh edges').opts(node_size=5)\n",
    "hnodes = circle.nodes.opts(size=5)\n",
    "dscirc = (hd.dynspread(hd.datashade(circle))*hnodes).relabel(\"Datashader edges\")\n",
    "\n",
    "circle + dscirc"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In these plots, you can select the Bokeh \"wheel zoom\" from the tool palette and then zoom using your scroll wheel or pan by clicking and dragging.  On a static site like a website or anaconda.org the lines in the Datashader plot will be rescaled blockily as you zoom in, but with a live server it will be dynamically re-rendered to show more detailed structure each time. \n",
    "\n",
    "You can try clicking and hovering on either plot to see what interactive features are available; in both cases the behavior for nodes should be the same (as the full set of nodes is being overlaid on both plots), while the edges also support interactivity in the pure-Bokeh version.\n",
    "\n",
    "As you can see, the pure-Bokeh version provides more interactivity, but the datashaded version will let you see the patterns of connectivity better for large graphs.  The datashader version will also work fine for  arbitrarily large graphs that would overwhelm the browser if used with Bokeh directly.  [HoloViews](http://holoviews.org/user_guide/Network_Graphs.html) makes it simple to switch between these two extremes as needed, using full-interactive plots for small datasets and adding whatever interactivity is required (as in the overlaid node plots on the right above) for larger datasets while rendering the full dataset as the main plot using datashader."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Real-world examples\n",
    "\n",
    "The above examples all use artificial datasets, but each real-world dataset will have its own specific properties. We've set up analyses of a few real datasets as well:\n",
    "\n",
    "### USA airport connections\n",
    "[![thumb](../assets/images/airport_connections.png)](https://anaconda.org/philippjfr/airport_connections)\n",
    "\n",
    "### Research institutions linked by joint UK grants\n",
    "\n",
    "[![thumb](../assets/images/uk_researchers.png)](https://anaconda.org/jbednar/uk_researchers)\n",
    "\n",
    "### PCAP computer network data\n",
    "\n",
    "[![thumb](../assets/images/pcap.png)](https://anaconda.org/philippjfr/packet_capture_graph_hv/notebook)"
   ]
  }
 ],
 "metadata": {
  "language_info": {
   "name": "python",
   "pygments_lexer": "ipython3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
